Gravitational instability of the primordial plasma: anisotropic evolution of structure seeds 
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Abstract 

We study how the presence of a background magnetic field, of intensity compatible with current observation constraints, affects 

^rj the linear evolution of cosmological density perturbations at scales below the Hubble radius. The magnetic field provides an 

additional pressure that can prevent the growth of a given perturbation; however, the magnetic pressure is confined only to the plane 

^ orthogonal the field. As a result, the "Jeans length" of the system not only depends on the wavelength of the perturbation but also 

*~i on its direction, and the perturbation evolution is anisotropic. We derive this result analytically and back it up with direct numerical 

00 integration of the relevant ideal magnetohydrodynamics equations during the matter-dominated era. Before recombination, the 

t-H kinetic pressure dominates and the perturbations evolve in the standard way, whereas after that time the magnetic pressure gets 

much larger and we observe the anisotropic evolution. We quantify this effect by estimating the eccentricity e of a Gaussian 

/"■n perturbation in coordinate space that was spherically symmetric at recombination. For a galactic-sized perturbation, we find that 



u 



e = 0.7 at z - 10 taking the background magnetic field of order 10 Gauss. 
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1. General Remarks 

Our theoretical knowledge of the Universe is based on 
the Standard Cosmological Model, that provides a convenient 
framework to satisfactorily explain the majority of cosmologi- 
cal observations, like the anisotropy pattern of the cosmic mi- 
crowave background (CMB) [1, 2, 3], the large-scale structure 
of the Universe [4, 5, 6], the Hubble diagram of distant type la 
supernovae [7, 8, 9] and the abundances of light elements [10]. 

The Standard Cosmological Model relies on the assumption 
that the Universe, at least at large scales, is highly homoge- 
neous and isotropic, and its geometry is thus described by the 
Robertson- Walker metric. In fact, the distribution of luminous 
red galaxies shows that the present Universe is homogeneous on 
scales greater than ~ 100 Mpc [11], while the isotropy of the 
CMB itself, which has a black-body distribution at T — 2.73 K 
with temperature fluctuations of order 10~ 5 or less, is an indi- 
cation of the isotropy of the Universe as a whole and a strong 
evidence for homogeneity at the time of hydrogen recombina- 
tion (nearly 400.000 years after the Big Bang, corresponding 
to a cosmological redshift z rec = 1100). On the other hand, 
below the "homogeneity scale" of 100 Mpc, the distribution of 
matter is definitely inhomogeneous. Such a dichotomy between 
the smoothness in the matter-energy distribution at z = Zrec an d 
the dumpiness of the recent Universe, z < 1 below a certain 
scale is explained by the mechanism of gravitational instabil- 
ity: the structures we observe today have been formed through 
the growth of tiny density perturbation seeds, that according to 
the currently accepted model were created in the early Universe 
during a phase of inflationary expansion. 



The presence of a large scale (i.e., coherent over a Hub- 
ble length), strong magnetic field is forbidden by the observed 
isotropy, as it would naturally single out a preferred spatial di- 
rection. However, a background, uniform magnetic field could 
be present at cosmological scales provided that its intensity is 
small enough. In particular, upper limits of order ~ 10~ 9 G 
on the present field intensity have been derived from obser- 
vations of the CMB temperature anisotropies [12, 13] and of 
its temperature-polarization correlation [14, l] 1 . Smaller scale 
fields, on the other hand, could be as strong as 10~ 6 G. 

The effects of large-scale magnetic fields on the evolution 
of cosmological structures have been studied extensively in the 
literature (for a complete review, see Ref. [16] and references 
therein), where both Newtonian and general relativistic treat- 
ments, the latter often using covariant and gauge-invariant tech- 
niques, are present. It is known that, among others, magnetic 
fields slow down (and possibly prevent) the growth of pertur- 
bations and can produce vorticities and shape distortions in the 
density field [17, 18]. In this paper, our goal is to revisit the is- 
sue of the existence of a "magnetic Jeans length" and of its de- 
pendence on the direction along which the perturbation propa- 
gates, other than to give a realistic estimate of its value in a dark 
matter-dominated Universe. The presence of a magnetic Jeans 
length has been discussed in several papers [16, 17, 18, 19, 20] 
both in a Newtonian and in general relativistic framework (for 
an analysis of the standard Jeans mechanism in presence of dis- 



It has been argued (see e.g., [15]) that the limits obtained from the CMB 
temperature data can be significantly relaxed in the presence of free-streaming 
neutrinos. However this does not affect the limits from the polarization. 



Preprint submitted to Physics Letter B 



July 19, 2011 



sipative effect see also [21, 22, 23]), but some of the analyses 
failed to recognize its angular dependence. Here, we present 
a neat derivation of the relevant instability scales bases on the 
equations of magnetohydrodynamics (MHD) on an expanding 
Universe. We discuss a simple generalization of the magnetic 
Jeans length suited for two-fluid systems, and clarify some mis- 
understandings that are present in the literature. We also show 
the results of the numerical integration of the coupled MHD and 
Poisson equations. Finally, we numerically study the distortion 
introduced in the density field by the anisotropy in the critical 
length. 

The paper is organized as follows. In Sec. 2, we characterize 
the Universe as a plasma. In Sec. 3, we introduce the basic 
equations, and in Sec. 4 we carry out the linearization pro- 
cedure. We derive the existence of the magnetic Jeans length 
from analytical considerations in Sec. 5, while in Sec. 6 we 
show some numerical results. Finally, we draw our conclusions 
in Sec. 7. 



2. The plasma features of the pre-recombination Universe 

Between the time of e + e~ annihilation (i.e., for a temperature 
T = 7> e - =* lMeV =* 2 x 10 10 K) and the time of hydro- 
gen recombination, roughly corresponding to the last scattering 
of the CMB (T = T KC =* 0.25 eV =* 3000 K), the cosmolog- 
ical fluid can be modeled as a plasma of electrons, positrons 
and photons. We do not consider here the neutrino component 
because for T < 1 MeV it evolves independently of the other 
components, and interacts only through gravity. The cosmo- 
logical baryon-to-photon ratio is extremely small and equal to 
ftb/fty - 6 X 10~ 10 , where n y and «/, are the photon and baryon 
number densities, respectively. Both ri/, and n y evolve with the 
redshift z according to the relation n^ 7 (z) = «2 (1 + z) 3 , where 
n° b =* 2.5 x 10 cm -3 and n y =* 410 cm -3 are their present val- 
ues. In general, in the following we will use a superscript 
to denote the present value of a quantity. Most of the baryons 
in the Universe are in the form of hydrogen nuclei, so that for 
simplicity we assume n& = n p (n p being the proton number 
density). Furthermore, n p is also equal to the electron number 
density n e - n p , because of the charge neutrality of the Uni- 
verse. 

A fundamental quantity characterizing a plasma is the Debye 
length Ad, namely the length over which electrons screen out 
electric fields in a plasma. It defines the length scale over which 
a system can consistently considered to be a plasma. The Debye 
length of a hydrogen plasma at temperature T is given by the 
relation Ad = ~\jT/4-nn e e 2 (e is the proton charge) and then, 
using the values given above: 



7l' = 2.73 K, and then 
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Before recombination, photons, protons and electrons share a 
common temperature T - T y that, in the absence of entropy- 
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Such a specific dependence on the redshift implies that the co- 
moving Debye length Ad = Ad(1 + z) is constant during the 
cosmological evolution and equal to Ad — 2 x 10 4 cm. 

Plasma effects can be relevant in a system when its physical 
dimension L is much larger than the Debye length. For the Uni- 
verse, the relevant length is the Hubble radius L = In = cH' 1 , 
where H is the Hubble parameter. This length represents the 
maximum scale at which microphysical processes can operate 
in order to establish the thermodynamical equilibrium. Today, 
Ih — 10 28 cm; during the matter dominated era l H oc (I + z)~ 3 ^ 2 , 
while in the radiation dominated era l H cc (1 + z)~ 2 . From the 
analysis of both these scales, it is evident that Ih » Ad turns 
out in the period considered here. This means that the cosmo- 
logical fluid can be considered as neutral at all scales of cosmo- 
logical interest. We also note that, under the same hypotheses, 
the baryonic matter Md within a Debye sphere is also constant 
and given by 



M D - | nm p nt,A D - 10 
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where m p is the proton mass. This "Debye mass" clearly results 
to be much smaller than any other of cosmological interest. 

Another meaningful index is the so-called plasma parameter 
Nd, i.e., the number of particles within a Debye sphere: 



N D = § nn h A D ■ 



(4) 



The dependence of Ad and n p on the redshift z implies that also 
Nd is a constant; in particular, since No - 10 7 » 1, the cosmo- 
logical fluid results to be a weakly coupled plasma. 

In order to provide a complete characterization of the cos- 
mological plasma, we now turn our attention to the plasma dis- 
sipative properties. Let us start by analyzing the redshift de- 
pendence of the plasma resistivity rj. For an electron-proton 
plasma, this is given by rj - m e v e j / n e e 2 , where m e is the elec- 
tron mass and v e i is the electron-ion-collision frequency. For 
the case under consideration, v e , is well approximated by the 
electron-electron collision frequency v ee [24], i.e., 



v,~ v „„ ~ 2.91 xl0~ 6 sec" 1 



lcm- 3 /leV/ 
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where In Ac is the Coulomb logarithm, introduced to quan- 
tify the effects that small-angle-diffusion collisions have in the 
Coulomb scattering. A simply estimate of Ac in a plasma is 
given by Ac - \2jiNd, so that for the cosmological fluid, the 
Coulomb logarithm is =* 20. Substituting Eq.(5) into the ex- 
pression for the resistivity given above, we get 



rj(z) -1.6x 



1+z 
1100 
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Close to recombination, i.e., z - Zdec — 1 100, the cosmological 



generating processes, evolves like T y (z) - T y (l + z), where plasma has an electric resistivity equal to Jj\ z=Zdec - 1.6 Ohm cm, 



i.e., a conductivity =* 0.6 Siemens cm ', a value typical of 
a semiconductor. In Gaussian units, used in the following, 



n\z 



1.8x10" 



Let us now focus on the viscous properties of the plasma. 
The shear viscosity coefficient of matter strongly coupled with 
radiation can be expressed as [25, 26] 
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a SB T t . 



(7) 



where a SB =* 5.7 x 10 _8 WK~ 4 rrr 2 is the Stefan-Boltzmann 
constant, while r denotes the mean collision time between par- 
ticles and can be estimated as t =* (n y cr T vy l (here, we have 
introduced cr T - 6.6 x 10~ 29 m 2 as the cross section for Thom- 
son scattering, and v =* c). For a photon gas at equilibrium at 
temperature T, n y =* 2.0 x 10 7 (r/K) 3 irr 3 and thus we finally 
obtain 

'T\ kg 



7] V (T)* 1.2x10" 



\K/ms 



(8) 



The resistivity and viscosity coefficients enter the MHD 
equations through the diffusion coefficients ij = rjc 2 /4n and 
fj v = rjy/po, where po i s assumed to be the background den- 
sity of the fluid, taken equal to the cosmological baryon den- 
sity, pb - 4.2 x 10~ 28 (1 + z) 3 kg m -3 (corresponding to a 
density parameter Q#/z 2 = 0.022). Using this and the relation 
T = r°(l + z), we get 
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The relative magnitude of the viscous and magnetic diffusion 
rates can be parameterized through the magnetic Prandtl num- 
ber Pr„ = ij v /fj. Using the expression above for i) v and fj we 

get 



Pr m ^5.0x10' 



1+z 

1 + Zrec 
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so that Pr m » 1, i.e., viscous diffusion is more important than 
resistive diffusion, at recombination and indeed always at the 
redshifts under consideration. 

Having discussed the relative importance of viscosity- and 
resistivity-driven dissipative effects, we now discuss at which 
scales these effects are relevant. In a magnetized plasma, a 
useful parameter is the Lundquist number S = Lv a /T], where 
L is a typical length scale and va is the Alfven velocity, re- 
lated to the ambient magnetic field Bo by v 2 A = B^/Anpa. The 
Lundquist number is basically the ratio between the resistive 
diffusion timescale t = L 2 /i) and the Alfven crossing timescale 
Ta - L/v a . Assuming for the magnetic field at recombina- 
tion the value Soiree - 10~ 3 G (this corresponds to a present 
intensity of 10~ 9 G) and, considering again po = ps, we obtain 
v A =* 1.2xl0 5 m/s[(l+z)/(l+z ree )] 1/2 . We thus find 



S =2.8x10 



1+z 



'(- M 
\Mpc/\l +z. 



(12) 



We also compare ta with the viscous diffusion timescale de- 
fined as r v = L 2 /fj v . The ratio 5 V = T v /ta can be thought as a 
viscous analogous to the Lundquist number: 



S v = 5.7xl0 9 



L 



1+z 



Mpc/\1 +Zrec 



(13) 



The baryonic mass (at the average background density) con- 
tained in a sphere of radius equal to the length where S ~ 1 can 
be calculated to be ~ 1O~ 51 M [(1 + z)/(l + z rec )] -3 , while the 
analogous quantity for S v is ~ 1O _1O M [(1 + z)/(l + z ree )] _9/2 . 
Both mass scales are well below the values of cosmological rel- 
evance at the redshifts of interest. 

The above considerations can be extended also to redshifts 
z < 1100. After that time most of the protons and electrons 
in the Universe exist inside the neutral atoms of light elements; 
however, a small residual ionized fraction x e = 2.5 x 10~ 4 sur- 
vives. The calculations carried on so far can be repeated for 
this post-recombination plasma, [taking also into account that 
the baryon temperature scales like (1 + z) 2 ] and the results are 
found to be qualitatively similar. 

Concluding, the analyses above show how the following hi- 
erarchy among the relevant time scales holds for all mass range 
and redshifts of interest: 



T A « T v «: T r 



(14) 



meaning that viscosity always dominates over resistivity, and 
that both dissipative effects can indeed be neglected when 
studying the propagation of Alfven waves. In view of the es- 
timation above, in the following we will regard the cosmolog- 
ical medium as an ideal plasma, implying that its dynamical 
description will be based on an ideal MHD picture. 

3. Basic equations 

In this section, we derive the basic equations describing the 
linear evolution of instabilities in the cosmological fluid, mod- 
eled as a magnetized plasma. In this respect, we underline that 
the full investigation of the perturbative dynamics of the Uni- 
verse would require a general-relativistic treatment, in order to 
correlate the matter and geometrical fluctuations. However, as 
long as one is interested in scales much smaller than the Hub- 
ble radius, i.e., L <K cH , a Newtonian treatment provides 
a consistent description of the dynamics. Nonetheless, in this 
scenario the expansion of the Universe can be accounted as the 
bulk background motion of the fluid [25, 27, 28]. 

The starting point of our treatment is the Eulerian set of 
equations governing the fluid motion, on which one can de- 
velop a perturbative theory by adding small fluctuations to the 
unperturbed cosmological background solution. The zeroth- 
order dynamics is derived by considering a flat homogeneous 
and isotropic Universe whose energy density is dominated by 
non-relativistic matter, and correctly describes the expansion of 
the Universe. We assume that a background magnetic field is 
present, whose contribution to the total energy density of the 
Universe can be considered negligible. 



Let us now start by briefly recalling the basic equations of 
non-relativistic, ideal and single fluid MHD, which govern the 
plasma motion. The mass conservation and the Newtonian 
gravitational field are described by the continuity and Poisson 
equations, the single-fluid dynamics is described by the Euler 
equation in presence of a magnetic field B and, finally, the elec- 
tromagnetic interaction can be summarized by the frozen-in and 
the Gauss laws. Such equations read 



d,p + V • pv - , 


(15a) 


V 2 (D - AnGp = , 


(15b) 


,v + p(v ■ V)v + VP+ 




+pVd> -(VxB)x B/4n = , 


(15c) 


a,B-Vx(vxB) = 0, 


(15d) 


VB = 0, 


(15e) 



respectively, where p is the mass density, v is the velocity field, 
<t> is the gravitational potential and G is Newton constant. This 
system constitutes the base of our perturbative approach. 

To derive the zeroth-order dynamics we assume the 
Robertson-Walker metric, i.e., ds 2 = dt 2 - a 2 (i)dl 2 , where 
a = a(t) represents the cosmological scale factor, and a per- 
fect fluid energy-momentum tensor as the matter source of the 
gravitational field, i.e., T M V = diag[po, -Po> _ Po, -Po], with 
Po = Po(0- In this scheme, the behavior of the mass density with 
time is obtained from the energy-momentum conservation-law 
TZ = and from the Friedmann equation, i.e., 



po + 3H(po + Po) = , 
a 2 +<K- §7rGpofl 2 = , 



(16) 

(17) 



respectively [the dot ( ' ) denotes the total derivative with respect 
to synchronous time]. Here H - a/a is the Hubbe parameter 
and 'K - const, is the curvature factor. 

Setting the matter-dominated Universe equation of state 
(EoS) Po ~ (Po <K po) in Eq.(16), the zeroth-order solution 
of the system (15) turns out to be 



P 

Po = -3 



v = Hr , B = — 



V<D 



f 7rGp r . 



(18) 



where p and Bo are dimensional constants, r (r — \r\) de- 
notes the radial coordinate vector and, of course, a(t) satisfies 
Eq.(17). We observe how this non-stationary solution char- 
acterizing the background dynamics is not affected by the so- 
called "Jeans swindle" proper of the static solution [25, 27, 28]. 
To obtain now the explicit time dependence of the unper- 
turbed quantities involved in the model, we restrict the analysis 
to the flat case, i.e., 'K = 0. From the Friedmann equation (17) 
and using the solution for po, one readily obtains 



= (6ttGp) 
1 



1/3 



2/3 



Po 



(19a) 
(19b) 



6nGt 2 ' 

Finally, we recall that the adiabatic sound speed is defined by 
v s = ^[dPjdp. For a general specific heat ratio y, we assume 



that the pressure varies as P = Kp y , so that the speed of sound 
is given by 



v^yKpV 



7 K r 2y+2 



(6nG)y- 1 



(19c) 



4. Perturbation scheme 



In order to analyze the implications that the physics of an 
ideal magnetized plasma can have on the structure formation, 
we will follow the standard perturbation approach. In this re- 
spect, we consider small perturbations around the zeroth-order 
cosmological solution derived above, i.e., we write p = po + p\ 
(with p\ <k po) and similarly for the other quantities P, v, <5 
and B. Substituting the perturbed quantities in Eqs.(15) and 
keeping only terms up to first order, one gets 



d,pi + 3Hp x + H(r ■ V)pi +p V • V\ = 
V 2 Oi - 4nGpi = 
d,v\ + Hv\ + H(r ■ V)vi + v 2 Vp 1 /p + 

+V«Di - (V x Bi) x B /(4;rpo) = 
d,Bi +2HBi +H(r-V)B l + 

+B (V • vi) - (Bo • V)vj = 

V-«! =0 



(20a) 
(20b) 

(20c) 

(20d) 

(20e) 



where, as already discussed, the pressure and density perturba- 
tions have been related through the adiabatic sound speed, i.e., 
P\ - v 2 p\. We are assuming that B^/4npo — v 2 A <K 1, where 
Bo = l-Bol> i n order to preserve the isotropy of the background 
flow. 

In the following, we replace B\ with the dimensionless mag- 
netic fluctuation b 1 = B\/Bq. Moreover, the analysis of the sys- 
tem above can be simplified by Fourier-transforming the spatial 
dependence of the involved quantities, i.e., using perturbations 
in the form of plane waves, taking 



Y {r,t) = Ht)e ik 



(21) 



where <pi - {p\, \>\, 4>i, B\) is the generic perturbation vari- 
able, and k is the physical wavenumber scaling as l/a(f). It is 
convenient to consider also the comoving wavenumber q - ak, 
that stays constant during the expansion. The evolution for 
a given harmonic can be simply obtained by the equations 
in real space with the substitutions <f>\ — > 0i, V — > ik and 
d t — » d t - iH(k ■ r). In the following, for the sake of simplicity, 
we will drop the tilde over the Fourier transformed variables. 
Then, the system (20) reduces to (hats denote unit vectors): 



i>i + Hv\ + i 



Po 



4nG 
It 



pi + 3Hpi + ip (k ■ vi) = , (22a) 



p x k + iv\ B x(kxbi) = 0, (22b) 



bi + iB Q (k ■ vi) - /(Bo • k)vt = , (22c) 



where we have used the Alfven velocity v 2 A - Bg/(47rpo), and 
we have already eliminated Oi by means of the Poisson equa- 
tion in A:-space, i.e., k 2 Q>\ - —4nGp\. It is understood that the 
constraint k ■ b \ =0 always hold. 



Decomposing now V\ in its components v\ and v| parallel and 
orthogonal to the direction of q respectively, i.e., v\ = v\q + v^ 
(where v\ ■ q - 0), and introducing the following scalar vari- 
ables: 



8 =Pi/po , 



= i(k ■ vi) = ikv\ , (23a) 

v = ik(v\-B Q ), (23b) 



we finally get a further simplified system: 

6+0=0 , 

6 + 2H6 - cold -co 2 A b = 0, 

b + (\ - 1?)6 - fiv = , 

v + 2flv + /UW^B = , 



where we have defined 

H = Bo • q , 



C0 2 A = v A k2 



co Q = v s k 



(24a) 
(24b) 
(24c) 
(24d) 



■4nG Po , (25) 



and, of course, ^ \i ^ 1 . We stress that ai 2 , is not positive 
defined. 



5. Evolution of the density contrast and conditions for col- 
lapse 

The form (24) of the evolution equations has the advantage 
that it clearly expresses the relationship between the physical 
quantities involved, other that being very well suited for numer- 
ical integration. Some further analytical insight can however be 
gained by reducing it to an unique higher order equation for the 
variable 6(t). 

Considering the case of a matter-dominated Universe, and 
using the explicit time dependence of the quantities involved in 
the model in that case, i.e., Eqs.(19), with some algebra one can 
derive the following fourth-order differential equation for (5(f): 

9f 4 <5 (4) + 60f 3 c5 (3) + [9A 2 + 76 + 9Air 2v ] t 2 c5 (2) + 

+ [(12A 2 + 8) + 12A 1 (l-3v)r 2v ] t6 m + (26) 

+ [-6A21? + 3Ai(3A 2/ u 2 + 12v 2 - 2v) r 2v ] (5 = 0, 

where 5^ denotes the I th derivative of 5 with respect to time, 
and we have defined the following constants: 



v = y-4/3, Ai = v L s k L t 



2.2 f 2y-2/3 



A 2 



2,2 



COA 



(27) 



We recall that y is the specific-heat ratio (P ~ p y ) and that 
y ^ 4/3, i.e., v Js 0. 

The most general solution of Eq.(26) for 6 is found to be the 
superposition of four independent solutions <5, (/ = 1 , . . . , 4), 
given by: 

Si = Ai f 2 n[{a u , fla); ( b lu b 2i , b 3i ); ~^^ ] , (28) 



where the A,'s are arbitrary integration constants, 
P J r q[(ai,...,a p y,(b\,...,bqy, z] denotes the generalized hy- 
pergeometric function of argument z, and 



x x = (-1 + V A -)/6 
x 3 = (-1 + Va7)/6 



x 2 = (-1 - VA_)/6 
x 4 = (-1 - VK7)/6 . 

A ± = 13 - 18A 2 + 6 ^OK 2 - 2) 2 + 24yu 2 A 2 . 



(29a) 
(29b) 

(29c) 



The constant coefficients a and b depend, in general, on v, A 2 , 
/i and we report their complete expressions in Appendix A. 

We are now interested in discussing the asymptotic behavior 
of the hypergeometric functions in the limit of very small or 
very large argument, i.e., Ai/4v 2 f 2v »l or <k1. As in the non- 
magnetic case discussed in Ref. [25], we restrict the analysis of 
the asymptotic behavior to the range ^ v ^ 1/3, i.e., for the 
standard regime 4/3 ^ y ^ 5/3. 

From the asymptotic expansion of the f functions in the case 
of large argument, i.e., Ax/4v 2 t 2v » 1, the density contrast al- 
ways shows a damped oscillating behavior with time. In fact, in 
this regime corresponding to "large" wavenumbers k » 2v/v,f, 
it always exists at least one asymptotic solution proportional 
(apart from the oscillating factor) to positive power of the argu- 
ment Aj /4v 2 f 2v (i.e., to negative power in time). Such a leading 
term depends on the value of v. 

On the other hand, the asymptotic expansion of the solu- 
tions (28) in the limit Ai/4v 2 f 2v — > 0, corresponding to "small" 
wavenumbers k <K 2v/v s t, can be written as 



di-t*' +0(Ai/4-v z t Zv ). 



(30) 



In order for the gravitational collapse to occur, at least one of the 
modes has to be growing, i.e., x, > 0. It is fairly easy to show 
that xx, x 2 and x\ are always negative, whereas the sign of X3 
depends on y. and A 2 . In particular, when fi + 0, x 3 is negative 
irregardless of the value of A 2 , while when yi — 0, X3 is positive 
if A 2 < 2/3. This means that density modes with k <K 2v/v s t 
always grow outside the plane orthogonal to the magnetic field 
(/i + 0), while on that plane (ji - 0) they are kept stable if 
k 2 > 2/3v 2 f 2 (so that A 2 > 2/3 ), i.e., if the magnetic field is 
strong enough. 

Remembering that po = l/6nGt 2 , the condition on the 
wavenumber related to A] can be rewritten as: 



k^kj = 



24nGv 2 p 



(31) 



that is substantially the same as the usual Jeans condition for 
gravitational instability. In fact, in the non-magnetic case, (to 
which our analysis reduces for coa = 0) this is the only criterion 
that separates the growing and the decaying modes [25]. 

In a similar way, the condition related to A 2 can be rewritten 
as: 



k^k A = 



4nGp 167!- 2 Gp 2 



Bl 



(32) 



Summarizing, we find that the presence of a background 
magnetic field introduces an anisotropy in the stability crite- 
rion. While outside the plane orthogonal to Bo, the stability 
of the perturbations is dictated only by the standard Jeans con- 
dition k ^ kj, on that plane the unstable modes are those for 
which the conditions k < kj and k < k& both hold 2 . In other 
words, if k,A < kj (basically equivalent to va > v s ), there are 
Jeans-unstable modes (those in the window k& < k < kj) that, 
in the orthogonal plane, are stabilized by the magnetic pressure. 
The window of stable modes gets wider for larger values of the 
ambient magnetic field, as expected. We underline that these 
results are qualitatively the same as those obtained for a static 
and uniform background [29]. A similar analysis was carried 
on by the authors of Ref. [17] obtaining a similar results. How- 
ever, their derivation contained a mistake when separating the 
real and imaginary components of the evolution equations [25]. 
For this reason they find a second-order differential equation 
instead than the fourth-order one discussed here. 

6. Numerical Analysis 

In the previous section, we have gained an important insight 
on the effect of a background magnetic field on the evolution 
of density perturbations. We now show some results obtained 
through the direct numerical integration of the differential sys- 
tem (24). 

6.1. Preliminaries 

We will focus in the period of the cosmological evolution that 
goes from the onset of matter domination (z - 3000) to the time 
of reionization (z =* 10). We start from matter domination be- 
caues, before that time, the growth of density perturbation was 
slowed down and practically freezed by the rapid expansion of 
the Universe. In the matter-dominated Universe, a oc t 2 ^ and 
H - 2/3f. We can ignore the presence of a dark energy com- 
ponent since this is sub-dominant until very recent times. The 
time period that we consider can be divided into two distinct 
phases, i.e., before and after the recombination of hydrogen oc- 
curring at z — 1100. Before recombination, the baryons are 
completely ionized and they are tightly coupled to photons, so 
that the total pressure of the fluid is given by radiation pressure. 
After recombination, most of the protons and electrons are in 
the form of neutral hydrogen atoms, leaving a small ionized 
fraction x e — 2.5 x 10~ 4 . Moreover, photons have decoupled so 
that the baryon pressure is given just by their kinetic pressure, 
dropping down by several orders of magnitude with respect to 
its pre-recombination value. 

In view of this, we take the speed of sound of the cosmolog- 
ical fluid before recombination to be [25] 



1 k B T b o- 
3m p + k B T b cr 



(33) 



2 It is easy to verify how the physical meaning of the condition is that the 
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Figure 1: Top panel: Alfven (red solid line) and sound (blue dashed line) speed 
as functions of redshift. The discontinuity at z = 1 100 corresponds to the re- 
combination on neutral hydrogen. Bottom panel: Magnetic (red solid line) and 
standard (blue dashed line) Jeans mass as a function of redshift. The mass 
contained inside the Hubble radius (black dotted line) is also shown for com- 
parison. 



where T b = T y is the common photon-baryon temperature and 
cr = 4assT i /3nbkB — 1.5 X 10 9 is the specific entropy. After 
recombination, we take instead 



5 k B T b 



(34) 



timescale for gravitational collapse r c ~ L/v c ~ -iJL 3 jGM 
shorter than both the acoustic and Alfven timescales r, ~ L/v, 
t a ~ L/v A ~ 1/kvA- 



l/ifeVj and 



where it should be taken into account that after recombination 
T b differs from T 7 and scales like cr 2 instead than a . Finally, 
unless otherwise stated, in order to calculate the Alfven speed 

va — JBq/4-7tpo we take a magnetic field with a present day in- 
tensity of 10~ 9 Gauss, while the density of charged particles is 
equal to the total baryon density p b = n b m p before recombina- 
tion, and equal to the ionized density fraction x e n b m p after that 
time. A plot of the Alfven and sound speeds as funtions of red- 
shift is shown in Figure 1 . The sudden drop in the sound speed 
at recombination is due to the sharp decrease of baryon pressure 
after photon decoupling, while the equally sudden increase in 
the Alfven speed at the same time is due to the decrease in the 
density of charged particles. 

In a detailed model, the presence of different uncoupled com- 
ponents making up the matter content of the Universe should be 
taken into account. In fact, most of the matter (~ 80%) is in the 
form of cold dark matter (CDM), interacting with the baryon- 
photon fluid only through the gravitational force. Moreover, 
after recombination, recombined and ionized hydrogen atoms 
should be effectively treated as distinct fluids, that again are 
coupled only by gravitation. Thus, a proper treatment should 



rely at least on a two- and three-fluid description before and 
after recombination, respectively. Here, we will make some 
approximations in order to be able to use the single fluid de- 
scription and still be able to draw meaningful conclusions about 
the perturbations in the baryonic component. First of all, we 
ignore the perturbations in the CDM component, in the pre- 
recombination era, this is justified by the fact that the large ra- 
diation pressure prevents baryons to fall into the potential wells 
created by CDM; this is known to be true in the non-magnetic 
case but we expect it to hold also in the case under consider- 
ation since, as seen in the last section, the magnetic field only 
acts to increase stability. After recombination, our treatment 
will rigorously remain valid for some time, i.e., only as long as 
the baryon density fluctuations do not catch up with those in the 
CDM component. 

For what concerns the modeling of the gravitational inter- 
action between the neutral and ionized baryon components af- 
ter recombination, we proceed in the following way. It is rea- 
sonable to assume that the newly recombined baryons inherit 
the fluctuations they had prior to recombination, so that at re- 
combination dm - <5hii (we use standard astronomical notation 
for ions). The neutral component obeys its own continuity and 
Euler equations (the latter obviously lacking the Lorentz force 
term) that are coupled to the corresponding equations for the 
ions by the Poisson equation, whose perturbed version will read 
V 2 Oi = 4-nG(Spm + ^Phii)- Given that the reionization fraction 
is very small, po,m s> po,Hii and then, since 6hi - <5hii, it is 
immediate to see that (6pm + ^Phii) - <5pm - Po.hi^hii- In 
other words, the Poisson equation rewrites: 



V Oi =4ttGpo,hi£hii , 



(35) 



so that the fluctuations in the neutral component do not explic- 
itly appear and we have reverted to a one-fluid approximation, 
although with a different source term for the gravitational field. 
It is easy to see that the only place where this source appears in 
the final equations is in the expression for o>q. Thus we can con- 
clude that the gravitational effect of the neutral component can 
be taken into account by the simple prescription that the density 
that appears inside oj^ should always be taken to mean the total 
(i.e., neutral + ionized) density of baryons. We remark again 
that this approach will be valid only as long as dm = £hii, so 
that we expect it to cease to provide a proper description some 
time after recombination. 

Before showing the results of the numerical integration, we 
illustrate in the right panel of Figure 1 the evolution of the 
standard Jeans wavenumber (31) and of its "magnetic" coun- 
terpart (32). Consistently with the discussion above, we use 
the total baryon density in both expressions, excepted when 
computing the Alfven speed, where we have to use the den- 
sity of the ionized fraction (in other words, barring constants, 

k\ °c PO,HlPO,HIl/#o = X eP\m / B o)' In ° rdel " t0 toke 0Ut the 

change in kj and Ic^ due to the expansion, we follow the con- 
vention to express the results in terms of the mass contained 
inside the corresponding length scales l/kj^- In particular, we 
consider the total baryonic mass (irrespective of the ionization 
state), contained inside a sphere of radius 2n/kj^. It can be 



seen that the window of modes that are made stable by the mag- 
netic field, i.e., those between the red dashed and the black solid 
line, spans, right after recombination, ten orders of magnitude 
in mass. 

As noted in the introduction, the existence of a magnetic 
Jeans length has been studied previously and all expressions for 
the critical wavenumber agree, apart from numerical factors, 
with expression (32). However, we argue that the numerical 
estimates of this quantity appearing in the literature are some- 
what inaccurate since it is often taken po inside the formula 
to be equal to the present critical density p c =* 9 x 10~ 27 kg 
trT 3 . This yields at the present time a magnetic Jeans length 
A A ~ \jk A ~ 10 kpc for B = 10 -9 G and A a ~ 1 Mpc for 
Bo = 10~ 7 . 3 These values would be accurate if a) the Universe 
were baryon-dominated, and b) all the baryons were in ionized 
form, but this situation is never realized in practice. One should 
instead take into account the actual density of baryons and, af- 
ter recombination, the fact that only a small fraction is ionized. 
Following our previous discussion, this is realized by the sub- 
stitution Pq — > x e pQ HI , that yields Aa ~ 100(10 4 ) Mpc for 
Bo = 10 ~ 9 (10~ 7 ) G (the latter value being even larger than the 
present Hubble radius). This approximation holds for a lim- 
ited time after recombination; once the baryons start to fall into 
the potential wells created by CDM, the same reasoning leads 
to the substitution Pq — » po.miPm where p m =* 0.135p e is the 
total matter density. Inserting the numerical values one gets 
A A ~ 40 (4 x 10 3 ) Mpc for B = 10~ 9 (10~ 7 ) G. 

6.2. Results 

We now discuss the results of the direct numerical integra- 
tion of Eqs.(24). The initial conditions for the integration have 
been chosen using the fact that power-law solutions for 6 can 
be found in the limit t — » (because in that limit both co\ and 
cOq are proportional to r 2 ). There are four distinct solutions 
of this kind, but only one corresponds to a growing mode. We 
have matched the initial conditions to the asymptotic growing 
solution at the initial time of integration. The latter has been 
chosen so that all the modes of interest were outside the hori- 
zon at that time. Even if the initial time falls in what would be 
the radiation-dominated era, nevertheless we always consider a 
matter-dominated Universe. All results have been normalized 
to the initial value of the density contrast. 

In Figure 2, we show the evolution of the density contrast for 
three different wavenumbers k - (1.7, 0.36, 0.08) Mpc -1 (nor- 
malized at the present time), i.e., for the following baryonic 
masses M = (1.7 x 10 11 , 1.7 x 10 13 , 1.7 x 1O 15 )M . These 
masses roughly correspond to the galactic, cluster and super- 
cluster mass scales, respectively. For each mode, we show the 
evolution in both the direction parallel to the background mag- 
netic field (,u = 1) and orthogonally to that direction (,u = 0). 



3 The latter value is also often said to be of the order of the scale of a galaxy 
cluster, while in effect it is closer to the scale of a galaxy as it can be seen by 
the fact that the mass enclosed inside a sphere of 1 Mpc radius at the critical 
density is ~ 10"M Q . The reason why a galaxy is much smaller than 1 Mpc is 
that it has detached from the Hubble flow and undergone non-linear evolution, 
so that its density is much larger than the cosmological average [27]. 
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Figure 2: Evolution of the dimensionless density perturbation & with redshift. The three panels show three different mass scales roughly corresponding, from left to 
right, to the galaxy, cluster and supercluster scales. In each panel, we show the evolution of perturbations orthogonal (red solid line) and parallel (blue dashed line) 
to the background magnetic field. 
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Figure 3: Equal density contours of an initially spherically symmetric Gaussian perturbation at the galactic scale in the x-y plane, at different times. The contours 
correspond to (0.1, 0.2, . . . , 0.9) times the central density. The background magnetic field is directed along the y axis and has an intensity Bq = 10~ 9 G. 



The present-day strength of the field is Bq = 10~ 9 G. In all 
cases, the perturbation is initially growing but then starts to os- 
cillate once the Jeans mass (that is growing) becomes larger 
than the mass of the perturbation. This happens earlier for 
smaller scales. In this phase the magnetic pressure does not 
play any role, as the much larger radiation pressure is actually 
providing the force that prevents the collapse. In fact, there is no 
difference in the evolution parallel and orthogonal to the field, 
as the radiation pressure is isotropic. The situation changes dra- 
matically after recombination, when the baryon pressure drops 
and it is instead the magnetic pressure that opposes the growth. 
However the magnetic pressure can only act in the plane or- 
thogonal to the field. Thus, perturbations in the direction of the 
field can grow unhindered, while the perturbations that are or- 
thogonal can be stabilized. As it can be seen from Figure 2, this 
happens basically for all perturbations below the supercluster 
scale. At z = 10, the relative growth of parallel perturbations 
with respect to orthogonal ones is of order 100, 10 and 3 for the 
galactic, cluster and supercluster scales, respectively. 

The fact that, after recombination, the evolution of the den- 
sity contrast in the presence of a magnetic field is different 



in different directions leads to the reasonable expectation that 
some degree of anisotropy will be generated even in initially 
symmetric structures. In order to show this, we consider a 
Gaussian density fluctuation with standard deviation <x in co- 
ordinate space at recombination: 



S(X, tree) = S(X = 0, trec)e 



(36) 



where the x are comoving coordinates centered at the maxi- 
mum of the perturbation. After Fourier-transforming, we sepa- 
rately evolve the different harmonics in momentum space using 
Eqs.(24) and we finally transform back to obtain the perturba- 
tion in coordinate space at a later time. In Figure 3 we show, 
for a background magnetic field directed along the y axis and 
with a present intensity of 1 nGauss, the evolution of a per- 
turbation with <x = 0.5 kpc at recombination (so that the 3cr 
region encloses a mass M =* 1.5 x 1O U M at the mean bary- 
onic density, i.e., roughly the mass of a galaxy). In particular, 
we show equal density contours at z = 1000, 100 and 10. It is 
evident from the figure how the perturbation becomes progres- 
sively squeezed along the direction orthogonal to the magnetic 
field. 
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Figure 4: Left panel: Eccentricity e of a Gaussian perturbation at the galactic scale as a function of redshift. The background magnetic field increases by a factor 10 
with each curve, starting from Bo = 10~ 12 G (bottom curve) up to Bo = 10 ' G (top curve). The two upper curves are only distinguishable at large redshifts. Right 
panel: the same as the left panel, but for a perturbation at the cluster scale. 



In order to quantify the anisotropy in the perturbation, we 
consider the isodensity contour corresponding to half the value 
at the peak and calculate its eccentricity e = ~\Jl - b 2 /a 2 , where 
a and b are the lengths of the semi-major and semi-minor axis 
of the contour, respectively. In Figure 4, we show how the ec- 
centricity changes with redshift; for the parameters used above, 
we get e ^ 0.7 at z = 10. 

7. Conclusions 

We have studied the effect of a background magnetic field 
on the linear evolution of cosmological density perturbations at 
scales well below the Hubble length, where a Newtonian treat- 
ment can be used, focusing on the matter-dominated era. The 
conditions that allow the growth of small density perturbations 
have been clearly stated; in particular, we have found that in the 
plane orthogonal to the ambient magnetic field, a new critical 
length appears, related to the presence of the magnetic pressure, 
while everywhere else outside that plane the stability is dictated 
by the standard Jeans criterion. This is also confirmed through 
a direct numerical integration of the relevant MHD equations 
during the matter-dominated era, where this effect is shown to 
be possibly important after recombination, when the magnetic 
pressure is much larger than the kinetic pressure of the resid- 
ual ionized nuclei. Finally, it has been shown how the depen- 
dence of the critical scale on the angle between the perturba- 
tion wavevector and the magnetic field could lead to a sizeable 
anisotropy in the perturbations at the onset of non-linearity. 

Our analysis has relied on some approximations: other than 
those already stated, we have ignored the gravitational effects 
of dark matter perturbation and we have assumed that after re- 
combination the perturbations in the neutral baryon component 
trace those in the charged component. We have argued that 
these approximations limit the validity of our treatment to some 



time after recombination. We defer a more detailed and fully 
general relativistic analysis, also taking into account the differ- 
ent fluid components, to a future work. 
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Appendix A. Hypergeometric Coefficients 

In following, we write the complete form of the coefficients 
of the hypergeometric function of Eq.(28). They read: 



a 1( i 2) = 1 + Va~ /12v - yJl-36v 2 A 2 /12v , 
ci(3 4) = 1 + Va7/12v - ^l-36yu 2 A 2 /12v , 
a 2( i 2) = 1 + Va~ /12v + yJl-36n 2 A 2 /12v , 
a 2( 3 4) = 1 + V^W12v + J 1 - 36yU 2 A 2 /12v , 



and 



bi(i 2 ) = 1 + VA_ /6v , 
fe 2 (i 2 ) = 1 + Va~ /12v - Va7/12v 
b 3( i 2 ) = 1 + a/A - /12v + VaT/12v 
fei (3 ) = 1 - a/a7/6v, 
2>2( 3 ) = 1- Va~ /12v- Va7/12v 
& 3(3 ) = 1+ Va~ /12v- Va7/12v 
ii ( 4) = l- Va~/12v+ Va7/12v 
b 2(i) = l+ Va~ /12v+ Va7/12v 
& 3(4 ) = 1 + VaT/6v. 
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